from scipy import log10,array,argsort,multiply,nan,inf,nanargmin,nanmax,isnan,concatenate,copy,shape
from pylab import imshow,xlabel,ylabel,xlim,plot,find
class FowlerNordheim(object):
  def __init__(self,V,I,dV=0.01):
    """ 
    Calculating Fowler Nordheim curves and minimum
    minimum are recorded if thay are dV from edge of voltage scan
    """
    self.x = copy(V)
    self.x[self.x==0] = nan
    self.x = 1./self.x
    self.y = log10(abs(multiply(self.x**2,I)))
    self.y[self.y==-inf]=nan
    self.x,self.y = self.x[self.x.argsort()],self.y[:,self.x.argsort()]
    self.x,self.y = [ self.x[self.x>0], -self.x[self.x<0] ], [ self.y[:,self.x>0], self.y[:,self.x<0] ]
    self.FN_min = array([[nan,nan]]*len(I))
    for i in [0,1]:
      i0 = array(map(nanargmin,self.y[i]))
      notNan =find(~isnan(i0))
      if len(notNan)>0:
        i0 = i0[notNan].astype('i').astype('i')
        tmp1 = (self.x[i][i0] > 1./(array(map( lambda dat: nanmax(V[~isnan(dat)] if len(V[~isnan(dat)])>0 else [nan]) ,I))-dV)[notNan])*(1./self.x[i][i0] > dV)
        self.FN_min[:,i][notNan][tmp1] = (1./self.x[i])[i0][tmp1]*(-1)**i
    self.FN_min_true = ~isnan(self.FN_min).any(1)
  def plot(self,twoD=True,imsize=(100,100),**kw): 
    if twoD:
       imshow(concatenate(self.y[::-1],1).T,origin=kw.get('origin','lower'),aspect=kw.get('aspect','auto'),**kw)
       x=concatenate(self.x[::-1],1)
    else:
      plot(self.x[0],self.y[0][:].T,'.',label='positive')
      plot(self.x[1],self.y[1][:].T,'.',label='negative')
      xlim(0,10);ylim(-8,-3)
      xlabel(r"1/V [$V^{-1}$]"); ylabel(r"$\ln(I/V^2)$")
  